library(magrittr)
library(openxlsx)
library(data.table)
library(RMySQL)
library(dplyr)
library(ggplot2)
library(ggridges)
library(plotly)
library(DT)
library(htmlTable)
library(kableExtra)
library(tidyr)
library(lemon)
library(patchwork)
library(readr)
#-----------------------------------------------------------------------------------#
# summary setup
#-----------------------------------------------------------------------------------#

prod.effects <- dbGetQuery(con, "SELECT chemical.*,
                          study.study_id,
                          study.processed,
                          study.study_type,
                          study.study_year,
                          study.study_source,
                          study.species,
                          study.strain_group,
                          study.admin_route,
                          study.admin_method,
                          study.substance_purity,
                          study.dose_start,
                          study.dose_start_unit,
                          study.dose_end,
                          study.dose_end_unit,
                          endpoint.endpoint_category,
                          endpoint.endpoint_type,
                          endpoint.endpoint_target,
                          endpoint.endpoint_id,
                          tg_effect.life_stage,
                          effect.effect_id,
                          effect.effect_desc,
                          effect.cancer_related,
                          tg.sex,
                           tg.n,
                          tg.generation,
                          dose.dose_level,
                          dose.conc,
                          dose.conc_unit,
                          dose.vehicle,
                          dtg.dose_adjusted,
                          dtg.dose_adjusted_unit,
                          dtg_effect.*
                          FROM 
                          ((((((((prod_toxrefdb_2_0.chemical 
                          INNER JOIN prod_toxrefdb_2_0.study ON
                          chemical.chemical_id=study.chemical_id)
                          INNER JOIN prod_toxrefdb_2_0.dose ON dose.study_id=study.study_id)
                          INNER JOIN prod_toxrefdb_2_0.tg ON tg.study_id=study.study_id)
                          INNER JOIN prod_toxrefdb_2_0.dtg ON tg.tg_id=dtg.tg_id AND
                          dose.dose_id=dtg.dose_id)
                          INNER JOIN prod_toxrefdb_2_0.tg_effect ON tg.tg_id=tg_effect.tg_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_0.dtg_effect ON
                          tg_effect.tg_effect_id=dtg_effect.tg_effect_id AND
                          dtg.dtg_id=dtg_effect.dtg_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_0.effect ON
                          effect.effect_id=tg_effect.effect_id)
                          INNER JOIN prod_toxrefdb_2_0.endpoint ON
                          endpoint.endpoint_id=effect.endpoint_id) 
                          where study.processed=1 ;") %>% 
  data.table()

prod.wo.effects <- dbGetQuery(con, "SELECT chemical.*,
                          study.study_id,
                          study.processed,
                          study.study_type,
                          study.study_year,
                          study.study_source,
                          study.species,
                          study.strain_group,
                          study.admin_route,
                          study.admin_method,
                          study.substance_purity,
                          study.dose_start,
                          study.dose_start_unit,
                          study.dose_end,
                          study.dose_end_unit,
                          endpoint.endpoint_category,
                          endpoint.endpoint_type,
                          endpoint.endpoint_target,
                          endpoint.endpoint_id,
                          tg_effect.life_stage,
                          effect.effect_id,
                          effect.effect_desc,
                          effect.cancer_related,
                          tg.sex,
                           tg.n,
                          tg.generation,
                          dose.dose_level,
                          dose.conc,
                          dose.conc_unit,
                          dose.vehicle,
                          dtg.dose_adjusted,
                          dtg.dose_adjusted_unit,
                          dtg.mg_kg_day_value,
                          dtg_effect.*
                          FROM 
                          ((((((((prod_toxrefdb_2_0.chemical 
                          LEFT JOIN prod_toxrefdb_2_0.study ON
                          chemical.chemical_id=study.chemical_id)
                          LEFT JOIN prod_toxrefdb_2_0.dose ON dose.study_id=study.study_id)
                          LEFT JOIN prod_toxrefdb_2_0.tg ON tg.study_id=study.study_id)
                          LEFT JOIN prod_toxrefdb_2_0.dtg ON tg.tg_id=dtg.tg_id AND
                          dose.dose_id=dtg.dose_id)
                          LEFT JOIN prod_toxrefdb_2_0.tg_effect ON tg.tg_id=tg_effect.tg_id)
                          LEFT JOIN prod_toxrefdb_2_0.dtg_effect ON
                          tg_effect.tg_effect_id=dtg_effect.tg_effect_id AND
                          dtg.dtg_id=dtg_effect.dtg_id)
                          LEFT JOIN prod_toxrefdb_2_0.effect ON
                          effect.effect_id=tg_effect.effect_id)
                          LEFT JOIN prod_toxrefdb_2_0.endpoint ON
                          endpoint.endpoint_id=effect.endpoint_id)
                          where study.processed=1 ;") %>% 
  data.table()

sbox.effects <- dbGetQuery(con, "SELECT chemical.preferred_name,
                          study.chemical_id,
                          study.study_id,
                          study.processed,
                          study.study_type,
                          study.study_year,
                          study.study_source,
                          study.species,
                          study.strain,
                          study.strain_group,
                          study.admin_route,
                          study.admin_method,
                          study.substance_purity,
                          study.dose_start,
                          study.dose_start_unit,
                          study.dose_end,
                          study.dose_end_unit,
                          study.study_citation,
                          endpoint.endpoint_category,
                          endpoint.endpoint_type,
                          endpoint.endpoint_target,
                          endpoint.endpoint_id,
                          tg_effect.life_stage,
                          tg_effect.effect_desc_free,
                          tg_effect.no_quant_data_reported,
                          effect.effect_id,
                          effect.effect_desc,
                          effect.cancer_related,
                          tg.sex,
                          tg.generation,
                          tg.n,
                          tg.dose_duration, 
                          tg.dose_duration_unit,
                          dose.dose_level,
                          dose.conc,
                          dose.conc_unit,
                          dose.vehicle,
                          dose.dose_comment,
                          dtg.dose_adjusted,
                          dtg.dose_adjusted_unit,
                          dtg.mg_kg_day_value,
                          dtg_effect.*
                          FROM ((((((((prod_toxrefdb_2_1.chemical 
                          INNER JOIN prod_toxrefdb_2_1.study ON
                          chemical.chemical_id=study.chemical_id)
                          INNER JOIN prod_toxrefdb_2_1.dose ON dose.study_id=study.study_id)
                          INNER JOIN prod_toxrefdb_2_1.tg ON tg.study_id=study.study_id)
                          INNER JOIN prod_toxrefdb_2_1.dtg ON tg.tg_id=dtg.tg_id AND
                          dose.dose_id=dtg.dose_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.tg_effect ON tg.tg_id=tg_effect.tg_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.dtg_effect ON
                          tg_effect.tg_effect_id=dtg_effect.tg_effect_id AND
                          dtg.dtg_id=dtg_effect.dtg_id)
                          INNER JOIN prod_toxrefdb_2_1.effect ON effect.effect_id=tg_effect.effect_id)
                          INNER JOIN prod_toxrefdb_2_1.endpoint ON
                          endpoint.endpoint_id=effect.endpoint_id) 
                          where study.processed=1 ;") %>% 
  data.table()

sbox.wo.effects <- dbGetQuery(con, "SELECT chemical.preferred_name,
                          study.chemical_id,
                          study.study_id,
                          study.processed,
                          study.study_type,
                          study.study_year,
                          study.study_source,
                          study.species,
                          study.strain,
                          study.strain_group,
                          study.admin_route,
                          study.admin_method,
                          study.substance_purity,
                          study.dose_start,
                          study.dose_start_unit,
                          study.dose_end,
                          study.dose_end_unit,
                          study.study_citation,
                          endpoint.endpoint_category,
                          endpoint.endpoint_type,
                          endpoint.endpoint_target,
                          endpoint.endpoint_id,
                          tg_effect.life_stage,
                          tg_effect.effect_desc_free,
                          tg_effect.no_quant_data_reported,
                          effect.effect_id,
                          effect.effect_desc,
                          effect.cancer_related,
                          tg.sex,
                          tg.generation,
                          tg.n,
                          tg.dose_duration, 
                          tg.dose_duration_unit,
                          dose.dose_level,
                          dose.conc,
                          dose.conc_unit,
                          dose.vehicle,
                          dose.dose_comment,
                          dtg.dose_adjusted,
                          dtg.dose_adjusted_unit,
                          dtg_effect.*
                          FROM ((((((((prod_toxrefdb_2_1.chemical
                          LEFT OUTER JOIN prod_toxrefdb_2_1.study ON
                          chemical.chemical_id=study.chemical_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.dose ON dose.study_id=study.study_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.tg ON tg.study_id=study.study_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.dtg ON tg.tg_id=dtg.tg_id AND
                          dose.dose_id=dtg.dose_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.tg_effect ON tg.tg_id=tg_effect.tg_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.dtg_effect ON
                          tg_effect.tg_effect_id=dtg_effect.tg_effect_id AND
                          dtg.dtg_id=dtg_effect.dtg_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.effect ON
                          effect.effect_id=tg_effect.effect_id)
                          LEFT OUTER JOIN prod_toxrefdb_2_1.endpoint ON
                          endpoint.endpoint_id=effect.endpoint_id) 
                          where study.processed=1 ;")

study <- dbGetQuery(con, "SELECT * from prod_toxrefdb_2_1.study where processed=1;") %>% 
  data.table()

chemical <- dbGetQuery(con, "SELECT * from prod_toxrefdb_2_1.chemical;")

sbox.pod <- dbGetQuery(con, "SELECT * from prod_toxrefdb_2_1.pod;") %>%
  data.table()

#-----------------------------------------------------------------------------------#
# v2.1 POD extremes: lowest loael/lel, highest noael/nel
#-----------------------------------------------------------------------------------#

# study-level 
sbox.pod.L <- sbox.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
    filter(!is.na(mg_kg_day_value)) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(study_id, sex, effect_profile_id,pod_type) %>%
   slice(which.min(mg_kg_day_value))
sbox.pod.N <- sbox.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(study_id, sex, effect_profile_id,pod_type) %>%
   slice(which.max(mg_kg_day_value))
sbox.pod.study <-rbind(sbox.pod.N, sbox.pod.L)
sbox.pod.study1 <-sbox.pod.study[sbox.pod.study$effect_profile_id==1,]  
sbox.pod.study1 <- distinct(sbox.pod.study1)
sbox.pod.study2 <-sbox.pod.study[sbox.pod.study$effect_profile_id==2,]  
sbox.pod.study2 <- distinct(sbox.pod.study2)

#chemical-level
sbox.pod.L.chem <- sbox.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(chemical_id, sex, effect_profile_id,pod_type) %>%
   slice(which.min(mg_kg_day_value))
sbox.pod.N.chem <- sbox.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(chemical_id, sex, effect_profile_id,pod_type) %>%
   slice(which.max(mg_kg_day_value))

sbox.pod.chem <-rbind(sbox.pod.N.chem, sbox.pod.L.chem)
sbox.pod.chem1 <-sbox.pod.chem[sbox.pod.chem$effect_profile_id==1,]  
sbox.pod.chem1 <- distinct(sbox.pod.chem1)
sbox.pod.chem2 <-sbox.pod.chem[sbox.pod.chem$effect_profile_id==2,]  
sbox.pod.chem2 <- distinct(sbox.pod.chem2)

1 Overview

The Toxicity Reference Database (ToxRefDB) serves as a resource for structured animal toxicity data for many retrospective and predictive toxicology applications. ToxRefDB contains in vivo study data from over 5900 guideline or guideline-like human health relevant studies for over 1100 chemicals. For more comprehensive documentation, consult the ToxRefDB v2.1 User Guide or past publications, such as Watford et al., 2019. The purpose of this release note is to summarize and explore added data related to the ToxRefDB v2.1 minor update.

The study types covered in ToxRefDB include the following repeat dose study designs utilizing various administration routes (predominantly oral): chronic (CHR; 1-2 year exposures depending on species and study design) conducted predominantly in rats, mice, and dogs; subchronic (SUB; 90 day exposures) conducted predominantly in rats, mice, and dogs; subacute (SAC; 14-28 day exposures depending on the source and guideline) conducted predominantly in rats, mice, and dogs; prenatal developmental (DEV) conducted predominantly in rats and rabbits; multigeneration reproductive toxicity studies (MGR) conducted predominantly in rats; reproductive (REP) toxicity studies conducted largely in rats; developmental neurotoxicity (DNT) studies conducted predominantly in rats; and a small number of studies with designs characterized as acute (ACU), neurological (NEU), or “other” (OTH).

Many of the studies (over 3,000) come from registrant-submitted toxicity studies known as data evaluation records (DERs) from the U.S. EPA’s Office of Pesticide Programs (OPP). The majority of chemicals in the database are therefore pesticides. Since 2009, continued curation efforts have expanded ToxRefDB to include toxicity studies from ten additional sources, including the National Toxicology Program (NTP), peer-reviewed primary research articles (OpenLit), and pharmaceutical pre-clinical toxicity studies (Pfizer, Sanofi, GSK, Merck), among others (RIVM, PMRA, unpublished and unassigned sources).

ToxRefDB serves as a resource for study design, quantitative dose response, and endpoint testing status information given guideline specifications from the US Environmental Protection Agency (US EPA) and the National Toxicology Program (NTP) headquartered at the National Institute of Environmental Health Sciences. The legacy and current data curation workflow is described in more detail in later sections. A important component of ToxRefDB is its controlled vocabulary for studies and effects observed for enhanced data quality.

The first version of ToxRefDB (ToxRefDB 1.0) was initially released as a series of spreadsheets, which are still available on EPA’s FTP site and referenced in FigShare (https://doi.org/10.23645/epacomptox.6062545.v1). ToxRefDB underwent significant updates that are described in the recent publication (Watford et al., 2019) and was released as ToxRefDB v2.0. ToxRefDB v2.0 and associated summary files can be found here: https://doi.org/10.23645/epacomptox.6062545.v3.

ToxRefDB v2.1 is a minor update of ToxRefDB v2.0 to correct issues discovered with the compilation script that caused some extracted values to not import properly from AccessDB curation files, such as failure to import some effects. The .sql export of ToxRefDB v2.1 will be available for public download on Clowder, available here: https://clowder.edap-cluster.com/files/62e15e98e4b055edffbfb8c6?dataset=61147fefe4b0856fdc65639b&space=&folder=62c5cfebe4b01d27e3b2d851. Although the overall number of studies and chemicals remains unchanged, the v2.1 update includes additional data as previously curated studies with extracted dose treatment groups and effects are now fully accessible. This added data can improve the utility of ToxRefDB as a resource for curated legacy in vivo information by providing more complete information of the past animal studies conducted. Moving forward, an application-driven workflow with the Data Collection Tool (DCT) will be utilized to create a more sustainable process for loading curated information to a database and support a more regular release cycle.

The .sql export of ToxRefDB v2.1 and summary files is available for public download here: https://clowder.edap-cluster.com/files/62c5d1d1e4b01d27e3b2d86a?dataset=61147fefe4b0856fdc65639b&space=&folder=62c5cfebe4b01d27e3b2d851. In addition to the SQL downloads, ToxRefDB information is also summarized with calculated point-of-departure values at the chemical and study level for inclusion in the summary-level database, the Toxicity Value Database (ToxValDB), which is accessible via the CompTox Chemicals Dashboard: https://comptox.epa.gov/dashboard/chemical-lists/TOXREFDB2. ToxRefDB v2.1 values will be incorporated in the next ToxValDB release.

1.1 Summary of v2.0 to v2.1 Changes

The following table details a summary of differences between ToxRefDB v2.0 and v2.1. To create this table, both versions of ToxRefDB were filtered to present only data where a full curation with guideline profile observations was complete. This is achieved using a ‘processed’ flag set to 1 within the study table. The ToxRefDB v2.1 minor update has recovered thousands of extracted values that failed to import properly from the original AccessDB curation files. The v2.1 update includes added data, which improves the utility of ToxRefDB as a resource for curated legacy in vivo information. BMDExpress was not run to calculate BMD values for v2.1, therefore BMD tables were dropped from the v2.1 schema. The stored procedure for POD calculation was updated, as described in later sections.

unique_prod_effects<- unique(prod.effects[,c("dtg_id", "endpoint_category", "endpoint_type", "endpoint_target", "effect_desc" )]) 
unique_sbox_effects<- unique(sbox.effects[,c("dtg_id", "endpoint_category", "endpoint_type", "endpoint_target", "effect_desc" )]) 

critical_prod_effects <- prod.effects[prod.effects$critical_effect==1,]  
critical_prod_effects<- unique(critical_prod_effects[,c("dtg_id", "endpoint_category", "endpoint_type", "endpoint_target", "effect_desc" )])
critical_sbox_effects <- sbox.effects[sbox.effects$critical_effect==1,]  
critical_sbox_effects<- unique(critical_sbox_effects[,c("dtg_id", "endpoint_category", "endpoint_type", "endpoint_target", "effect_desc" )]) 


Output <- c("Total number of studies with complete curation", 
            "Number of studies with extracted effects",
            "Total number of chemicals", 
            "Total database rows, including studies with no extracted effects", 
            "Total effects extracted", 
            "Dose treatment groups with effects", 
            "Unique effects: Cholinesterase endpoint category", 
            "Unique effects: Developmental endpoint category", 
            "Unique effects: Reproductive endpoint category", 
            "Unique effects: Systemic endpoint category",
            "Unique critical effects: Cholinesterase endpoint category", 
            "Unique critical effects: Developmental endpoint category", 
            "Unique critical effects: Reproductive endpoint category", 
            "Unique critical effects: Systemic endpoint category")

v2.0 <- c(length(unique(prod.wo.effects$study_id)),
          length(unique(prod.effects$study_id)), 
          length(unique(prod.wo.effects$chemical_id)),
          nrow(prod.wo.effects), 
          nrow(prod.effects),
          length(unique(prod.effects$dtg_id)), 
          nrow(unique_prod_effects[endpoint_category=="cholinesterase",]), 
          nrow(unique_prod_effects[endpoint_category=="developmental",]),
          nrow(unique_prod_effects[endpoint_category=="reproductive",]),
          nrow(unique_prod_effects[endpoint_category=="systemic",]), 
          nrow(critical_prod_effects[endpoint_category=="cholinesterase",]), 
          nrow(critical_prod_effects[endpoint_category=="developmental",]),
          nrow(critical_prod_effects[endpoint_category=="reproductive",]),
          nrow(critical_prod_effects[endpoint_category=="systemic",]))

v2.1 <- c(length(unique(sbox.wo.effects$study_id)),
          length(unique(sbox.effects$study_id)), 
          length(unique(sbox.wo.effects$chemical_id)),
          nrow(sbox.wo.effects), 
          nrow(sbox.effects),  
          length(unique(sbox.effects$dtg_id)), 
          nrow(unique_sbox_effects[endpoint_category=="cholinesterase",]), 
          nrow(unique_sbox_effects[endpoint_category=="developmental",]),
          nrow(unique_sbox_effects[endpoint_category=="reproductive",]),
          nrow(unique_sbox_effects[endpoint_category=="systemic",]), 
          nrow(critical_sbox_effects[endpoint_category=="cholinesterase",]), 
          nrow(critical_sbox_effects[endpoint_category=="developmental",]),
          nrow(critical_sbox_effects[endpoint_category=="reproductive",]),
          nrow(critical_sbox_effects[endpoint_category=="systemic",]))

Change <- v2.1-v2.0
Table <- data.frame(Output, v2.0, v2.1, Change)

datatable(Table,
          filter='top', 
          options=list(pageLength = 15,searching=FALSE, autoWidth=FALSE,  scrollX=TRUE, initComplete = JS(
    "function(settings, json) {",
    "$('body').css({'font-family': 'Calibri'});",
    "}"
  )))

Treatment-related effects define the lowest effect levels (LELs) and no effect levels (NELs); these effects are simply treatment-related and significantly different from control. In contrast, toxicological opinion informs selection of a critical effect designation to define the lowest observable adverse effect and no observable adverse effect levels (LOAELs, NOAELs).

2 Data Landscape of v2.1

ToxRefDB v2.1 contains summary information from 5986 studies for 1143 chemicals. As part of the ToxRefDB 2.0 update, quantitative (i.e. dose-response) data was extracted. Currently, this curation has been completed for 3871 studies (indicated with processed=1) with plans to extract and release the remaining data in subsequent data releases.

To provide the reader with a summary of the scope and coverage of the database, ToxRefDB was filtered to present only data where a full curation with guideline profile observations was complete. This is achieved using a ‘processed’ flag set to 1 within the study table.

2.1 Summary Statistics

A summary table of the number of chemicals and number of studies for each study source, study type, and species

#-----------------------------------------------------------------------------------#
# Create summary statistics table of study information in toxrefdb v2.1 
#-----------------------------------------------------------------------------------#

#### select and group target information for table
sbox.study <- study
counts <- sbox.study 

#### adjust level names of study_source to reflect published table's names
# rename every study source which begins with "pharma_" to "Pharma"
counts$study_source<-gsub("pharma[^ ]*", "Pharma", counts$study_source)
# change other study_source names
'%notin%' <- Negate('%in%') 
counts <- counts %>% 
  mutate(study_source = replace(study_source, study_source == "ntp", "NTP")) %>%
  mutate(study_source = replace(study_source, study_source == "open_lit", "OpenLit")) %>%
  mutate(study_source = replace(study_source, study_source == "opp_der", "OPP DER")) %>%
  mutate(study_source = replace(study_source, study_source %notin% c("Pharma", "NTP", "OpenLit", "OPP DER"), "Other"))

counts <- counts %>%
  group_by(study_type, study_source, species) %>%
  summarize(n_studies = n_distinct(study_id), n_chems = n_distinct(chemical_id))
  

#### for each study type, extract its rows and calculate its totals for studies and chemicals

#subset counts by study type
CHR <- as.data.frame(counts[which(counts$study_type=="CHR"),])
DEV <- as.data.frame(counts[which(counts$study_type=="DEV"),])
MGR <- as.data.frame(counts[which(counts$study_type=="MGR"),])
SAC <- as.data.frame(counts[which(counts$study_type=="SAC"),])
SUB <- as.data.frame(counts[which(counts$study_type=="SUB"),])

# vector of study type names
s.type.char <- as.vector(unique(counts$study_type)) 
# counts data broken into a list by study type. Calling as.data.frame on an 
# element of this list will produce a data frame for all data from counts for the
# study type which corresponds to that element.
s.type.obj <- lapply(s.type.char, function(s) eval(parse(text=s)))  
# character vector of object names for the totals rows of each study type
totName<-paste(s.type.char, "t", sep=".")
# empty data frame to populate with data frames and totals rows of each study type
counts2<-as.data.frame(NULL)

for(i in 1:length(s.type.char)){
  # create a totals row for each study type with totals for studies and chemicals
  assign(totName[i], 
         c(paste("Total", s.type.char[i], sep = " "), 
           "", 
           "",
           n_distinct(sbox.study[which(study_type==s.type.char[i]),]$study_id), 
           n_distinct(sbox.study[which(study_type==s.type.char[i]),]$chemical_id)), 
         .GlobalEnv)
  # append each study's count data and totals row to a data frame combining those of all studies
  counts2 <- rbind(counts2, rbind(as.data.frame(s.type.obj[i]), eval(parse(text=totName[i])))) 
}

# create row for database totals for studies and chemicals, then append to counts2 data frame
dbtot <- c("Database totals", "", "", n_distinct(sbox.study$study_id), n_distinct(sbox.study$chemical_id))
counts2 <- rbind(counts2, dbtot)        

#### build table of combined data frame
table_dt <- counts2
kbl(table_dt, 
    col.names = c("Study type", "Study source", "Species", "Number of studies", "Number of chemicals"), 
    align = "c") %>%
  kable_paper(full_width = F) %>%
  column_spec(1, bold = T) %>%
  collapse_rows(columns = 1:2, valign = "top") %>%
  footnote(general = "CHR = Chronic, DEV = Prenatal-Developmental, MGR = Multigeneration Reproductive, SAC = Subacute, SUB = Subchronic", general_title = "")
Study type Study source Species Number of studies Number of chemicals
CHR NTP mouse 178 173
rat 169 164
OpenLit mouse 4 4
rat 5 5
OPP DER dog 331 298
hamster 4 3
mouse 342 303
primate 1 1
rat 398 328
Total CHR 1432 557
DEV NTP mouse 1 1
rabbit 3 3
rat 6 6
OpenLit rat 1 1
OPP DER mouse 18 16
rabbit 431 372
rat 508 433
Other mouse 1 1
rabbit 1 1
rat 4 4
Total DEV 974 486
MGR OpenLit rat 1 1
OPP DER mouse 2 2
rat 339 310
Other rat 19 19
Total MGR 361 331
SAC NTP mouse 29 26
rat 30 29
OPP DER dog 1 1
mouse 3 3
rabbit 6 6
rat 15 13
Total SAC 84 51
SUB NTP hamster 1 1
mouse 119 107
rat 127 114
OpenLit mouse 2 2
rat 4 4
OPP DER dog 214 195
hamster 4 4
mouse 123 112
primate 3 3
rabbit 5 4
rat 418 335
Total SUB 1020 498
Database totals 3871 748
CHR = Chronic, DEV = Prenatal-Developmental, MGR = Multigeneration Reproductive, SAC = Subacute, SUB = Subchronic

2.2 Study-Level Visuals

A breakdown of studies by study source, study type, and species

#### select and group target information for table
counts <- study 

#### adjust level names of study_source to reflect published table's names
# rename every study source which begins with "pharma_" to "Pharma"
counts$study_source<-gsub("pharma[^ ]*", "Pharma", counts$study_source)
# change other study_source names
'%notin%' <- Negate('%in%') 
counts <- counts %>% 
  mutate(study_source = replace(study_source, study_source == "ntp", "NTP")) %>%
  mutate(study_source = replace(study_source, study_source == "open_lit", "OpenLit")) %>%
  mutate(study_source = replace(study_source, study_source == "opp_der", "OPP DER")) %>%
  mutate(study_source = replace(study_source, study_source %notin% c("Pharma", "NTP", "OpenLit", "OPP DER"), "Other"))

#### dataframe with study type counts
counts.st <- counts %>%
  group_by(study_type) %>%
  summarize(n_studies = n_distinct(study_id), n_chems = n_distinct(chemical_id))

#### dataframe with species counts
counts.sp <- counts %>%
  group_by(species) %>%
  summarize(n_studies = n_distinct(study_id), n_chems = n_distinct(chemical_id))

#### dataframe with study source counts
counts.ss <- counts %>%
  group_by(study_source) %>%
  summarize(n_studies = n_distinct(study_id), n_chems = n_distinct(chemical_id))

#### study-level (study counts) frequency bar plot for study type, species, and study source

study.st <- ggplot(data=counts.st, aes(x=reorder(factor(study_type), -n_studies), y=n_studies)) +
  geom_bar(stat="identity", fill="aquamarine3") +
  theme_minimal() +
  labs(title="Type", x="study type", y="number of studies") +
  geom_text(aes(label=n_studies), vjust = -0.3, size=2) +
  ylim(0, max(counts.st$n_studies)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=0, vjust=1.2, hjust=1),
        aspect.ratio = 4/2)

study.sp <- ggplot(data=counts.sp, aes(x=reorder(factor(species), -n_studies), y=n_studies)) +
  geom_bar(stat="identity", fill="aquamarine3") +
  theme_minimal()  +
  labs(title="Species", x="species", y="number of studies") +
  geom_text(aes(label=n_studies), vjust = -0.3, size=2) +
  ylim(0, max(counts.sp$n_studies)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=45, hjust=1),
        aspect.ratio = 4/2)

study.ss <- ggplot(data=counts.ss, aes(x=reorder(factor(study_source), -n_studies), y=n_studies)) +
  geom_bar(stat="identity", fill="aquamarine3") +
  theme_minimal() +
  labs(title="Source", x="study source", y="number of studies") +
  geom_text(aes(label=n_studies), vjust = -0.3, size=2) +
  ylim(0, max(counts.ss$n_studies)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=45, hjust=1),
        aspect.ratio = 4/2)

#### all study-level plots together
study.st + study.sp + study.ss

2.3 Chemical-Level Visuals

A breakdown of chemicals by study source, study type, and species

#### chemical-level (chemical counts) frequency bar plot
chem.st <- ggplot(data=counts.st, aes(x=reorder(factor(study_type), -n_chems), y=n_chems)) +
  geom_bar(stat="identity", fill="darksalmon") +
   theme_minimal() +
  labs(title="Type", x="study type", y="number of chemicals") +
  geom_text(aes(label=n_chems), vjust = -0.3, size=2) +
  ylim(0, max(counts.st$n_chems)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=0, vjust=1.2, hjust=1),
        aspect.ratio = 4/2)

chem.sp <- ggplot(data=counts.sp, aes(x=reorder(factor(species), -n_chems), y=n_chems)) +
  geom_bar(stat="identity", fill="darksalmon") +
   theme_minimal() +
  labs(title="Species", x="species", y="number of chemicals") +
  geom_text(aes(label=n_chems), vjust = -0.3, size=2) +
  ylim(0, max(counts.sp$n_chems)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=45, hjust=1),
        aspect.ratio = 4/2)

chem.ss <- ggplot(data=counts.ss, aes(x=reorder(factor(study_source), -n_chems), y=n_chems)) +
  geom_bar(stat="identity", fill="darksalmon") +
   theme_minimal() +
  labs(title="Source", x="study source", y="number of chemicals") +
  geom_text(aes(label=n_chems), vjust = -0.3, size=2) +
  ylim(0, max(counts.ss$n_chems)*1.1) +
  theme(plot.title = element_text(size=10),
        axis.text.x=element_text(angle=45, hjust=1),
        aspect.ratio = 4/2)
 
#### all chemical-level plots together
chem.st + chem.sp + chem.ss

2.4 Chemical List

A list of chemicals included in ToxRefDB with associated study types

inspect <- distinct(study, chemical_id, study_type)
chemlist <- inspect %>% group_by(chemical_id) %>% summarise(val=paste(study_type, collapse=", "))
chemlist2 <- merge(chemical, chemlist, by="chemical_id")

chemlist2$chemical_id <- NULL
setnames(chemlist2, "val", "study_type")

datatable(chemlist2, 
          filter='top', 
          options=list(pagelength=25, autoWidth=FALSE,  scrollX=TRUE, initComplete = JS(
    "function(settings, json) {",
    "$('body').css({'font-family': 'Calibri'});",
    "}"
  )))

3 Reviewing Calculated Points of Departure Using Added Data

For purposes of predictive toxicology, points of departure (PODs) can be computed per chemical (i.e., lowest dose that produced effects or adverse effects across all study types included in the database) or per study (i.e., lowest dose that produced effects or adverse effects in a given study of interest).

Treatment-related effects define the lowest effect levels (LELs) and no effect levels (NELs); these effects are simply treatment-related and significantly different from control. In contrast, toxicological opinion informs selection of a critical effect designation to define the lowest observable adverse effect and no observable adverse effect levels (LOAELs, NOAELs).

In the ToxRefDB pod table, effects are grouped together within “effect profiles” for the purposes of POD derivation. The first effect profile calculates POD values for each study type, life stage, and endpoint category combination. A second effect profile calculates POD values for each endpoint category-endpoint type pairing, except in the case of the systemic endpoint category, where PODs were reported for each endpoint target (e.g., organs).

For these release note visuals, one set of “extreme” POD values (lowest loael/lel and highest noael/nel mg/kg/day value) are selected for each study id at the study-level, and for each chemical id at the chemical-level.

The stored procedure for POD calculation and resulting table has been updated for v2.1 with new columns (including sex, admin_route, and species) to support evolving needs and ToxValDB use.

3.1 Study-Level PODs

These visualizations show ToxRefDB v2.1 study-level PODs calculated across effect profile groups 1 and 2, plotted by sex.
95% of the effect profile #1 PODs lie in the gray rectangular region between 0.41 and 1000 mg/kg/day.
95% of the effect profile #2 PODs lie in the gray rectangular region between 0.7 and 1239.05 mg/kg/day.

study1 <- ggplot(sbox.pod.study1, 
                 aes(x = pod_value, 
                     y = reorder(study_id, pod_value), 
                     shape = pod_type, 
                      color = factor(pod_type, levels = c("lel", "loael", "nel", "noael")))) + 
  annotate("rect", 
           xmin = quantile(sbox.pod.study1$pod_value, 0.05), 
           xmax = quantile(sbox.pod.study1$pod_value, .95), 
           ymin = -Inf, 
           ymax = Inf, 
           alpha = 0.5, 
           fill = "lightgray") + 
  geom_point(position = position_jitter(h = 0.1),
             alpha = 0.5) + 
  theme_bw() +
  scale_shape_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c(17, 18, 15, 16)) +
  scale_color_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c("aquamarine3","darkviolet", "chartreuse", "darksalmon")) + 
  scale_x_log10(breaks = 10^pretty(log10(sbox.pod.study1$pod_value)), 
                labels = formatC(10^pretty(log10(sbox.pod.study1$pod_value)), 
                                 format = "e", 
                                 digits = 2)) +
  scale_y_discrete(position = "right") +
  labs(title = "Effect Profile #1", 
       x = "Conc (mg/kg/day)") +
  theme(panel.grid.major.x = element_line(color = "gray"),
        panel.grid.minor.x = element_line(color = "gray"),
        panel.grid.major.y = element_blank(),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom", 
        legend.box = "vertical",
        axis.text.y = element_blank(), 
        axis.title.y = element_text(size = 12),
        axis.ticks.y = element_blank(),
        axis.title.x = element_text(size = 12),
        legend.key.size = unit(1, 'cm'), #change legend key size
        legend.key.height = unit(1, 'cm'), #change legend key height
        legend.key.width = unit(1, 'cm'),
        legend.text = element_text(size = 12)) +
  facet_rep_wrap( ~ sex, 
                  ncol = 1, 
                  repeat.tick.labels = "bottom")




study2 <- ggplot(sbox.pod.study2, 
                 aes(x = pod_value, 
                     y = reorder(study_id, pod_value), 
                     shape = pod_type, 
                     color = factor(pod_type, levels = c("lel", "loael", "nel", "noael")))) +  
  annotate("rect", 
           xmin = quantile(sbox.pod.study2$pod_value, 0.05), 
           xmax = quantile(sbox.pod.study2$pod_value, .95), 
           ymin = -Inf, 
           ymax = Inf,
           alpha = 0.5, 
           fill = "lightgray") +
  geom_point(position = position_jitter(h = 0.1),
             alpha = 0.5) + 
  theme_bw()+
  scale_shape_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c(17, 18, 15, 16),
                     guide = FALSE) +
  scale_color_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c("aquamarine3","darkviolet", "chartreuse","darksalmon"),
                     guide = FALSE) +
  scale_x_log10(breaks = 10^pretty(log10(sbox.pod.study2$pod_value)), 
                labels = formatC(10^pretty(log10(sbox.pod.study2$pod_value)), 
                                 format = "e", 
                                 digits = 2)) +
  scale_y_discrete(position = "right") +
  labs(title = "Effect Profile #2", 
       x = "Conc (mg/kg/day)", 
       y = "Study" ) +
  theme(panel.grid.major.x = element_line(color = "gray"),
        panel.grid.minor.x = element_line(color = "gray"),
        panel.grid.major.y = element_blank(),
        panel.grid.minor.y = element_blank(),
        legend.position = "none", 
        legend.box = "none", 
        axis.title.y = element_text(size = 12),
        axis.ticks.y = element_blank(),
        axis.title.x = element_text(size = 12), 
        axis.text.y = element_blank())  +
  facet_rep_wrap( ~ sex, 
                  ncol = 1, 
                  repeat.tick.labels = "bottom") 
      

#### all plots together
study1 + study2 + plot_layout(guides = "collect") & 
  theme(legend.position = "bottom")  & 
  labs(x = "Conc (mg/kg/day)", 
       y = "Study")

3.2 Chemical-Level PODs

These visualizations show ToxRefDB v2.1 chemical-level PODs calculated across effect profile groups 1 and 2, plotted by sex.
95% of the effect profile #1 PODs lie in the gray rectangular region between 0.3 and 1027.5 mg/kg/day.
95% of the effect profile #2 PODs lie in the gray rectangular region between 0.3 and 1703 mg/kg/day.

chem1 <- ggplot(sbox.pod.chem1, 
                aes(x = pod_value, 
                    y = reorder(chemical_id, pod_value), 
                    shape = pod_type, 
                    color = factor(pod_type, levels = c("lel", "loael", "nel", "noael")))) + 
  annotate("rect", 
           xmin = quantile(sbox.pod.chem1$pod_value, 0.05), 
           xmax = quantile(sbox.pod.chem1$pod_value, .95), 
           ymin = -Inf, 
           ymax = Inf, 
           alpha = 0.5, 
           fill = "lightgray") +
  geom_point(position = position_jitter(h = 0.1),
             alpha = 0.5) + 
  theme_bw() +
  scale_shape_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c(17, 18, 15, 16)) +
  scale_color_manual("POD type", 
                     labels = c("LEL", "LOAEL", "NEL", "NOAEL"), 
                     values = c("aquamarine3","darkviolet", "chartreuse", "darksalmon")) + 
  scale_x_log10(breaks = 10^pretty(log10(sbox.pod.chem1$pod_value)), 
                labels = formatC(10^pretty(log10(sbox.pod.chem1$pod_value)), 
                                 format = "e", 
                                 digits = 2)) +
  scale_y_discrete(position = "right") +
  labs(title = "Effect Profile #1", 
       x = "Conc (mg/kg/day)", 
       y = "Chemical") +
  theme(panel.grid.major.x = element_line(color = "gray"),
        panel.grid.minor.x = element_line(color = "gray"),
        panel.grid.major.y = element_blank(),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom", 
        legend.box = "vertical", 
        axis.title.y = element_text(size = 12), 
        axis.ticks.y = element_blank(),
        axis.title.x = element_text(size = 12), 
        axis.text.y = element_blank(),
        legend.key.size = unit(1, 'cm'), #change legend key size
        legend.key.height = unit(1, 'cm'), #change legend key height
        legend.key.width = unit(1, 'cm'),
        legend.text = element_text(size = 12)) +
  facet_rep_wrap( ~ sex,
                  ncol = 1, 
                  repeat.tick.labels="bottom") 




chem2 <- ggplot(sbox.pod.chem2, 
                aes(x = pod_value, 
                    y = reorder(chemical_id, pod_value), 
                    shape = pod_type, 
                    color = factor(pod_type, levels = c("lel", "loael", "nel", "noael")))) +
  annotate("rect", 
           xmin = quantile(sbox.pod.chem2$pod_value, 0.05), 
           xmax = quantile(sbox.pod.chem2$pod_value, .95), 
           ymin = -Inf, 
           ymax = Inf, 
           alpha = 0.5, 
           fill = "lightgray") +
  geom_point(position=position_jitter(h = 0.1),
             alpha = 0.5) + 
  theme_bw() +
  scale_shape_manual(values = c(17, 18, 15, 16), 
                     guide = FALSE) +
  scale_color_manual(values = c("aquamarine3","darkviolet", "chartreuse",  "darksalmon"), 
                     guide = FALSE) +
  scale_x_log10(breaks = 10^pretty(log10(sbox.pod.chem2$pod_value)), 
                labels = formatC(10^pretty(log10(sbox.pod.chem2$pod_value)), 
                                 format = "e", 
                                 digits = 2)) +
  scale_y_discrete(position = "right") +
  labs(title = "Effect Profile #2", 
       x = "Conc (mg/kg/day)", 
       y = "Chemical" ) +
  theme(panel.grid.major.x = element_line(color = "gray"),
        panel.grid.minor.x = element_line(color = "gray"),
        panel.grid.major.y = element_blank(),
        panel.grid.minor.y = element_blank(),
        legend.position = "none", 
        legend.box = "none", 
        axis.title.y = element_text(size = 12),
        axis.ticks.y = element_blank(),
        axis.title.x = element_text(size = 12), 
        axis.text.y = element_blank())  +
  facet_rep_wrap( ~ sex, 
                  ncol = 1, 
                  repeat.tick.labels = "bottom") 
      
#### all plots together
chem1 + chem2 + plot_layout(guides = "collect") & 
  theme(legend.position = "bottom")  & 
  labs(x = "Conc (mg/kg/day)", 
       y = "Chemical")

3.3 Study-Level Changes in v2.0 to v2.1 PODs

There are 3614 studies with calculated POD values in v2.1 to compare back to v2.0. The v2.0 pod table included calculated POD values for 5922 studies: 5904 with 4 POD types and 19 with less than 4 POD types. To do the comparison, the v2.1 POD calculation was rerun against v2.0 schema since calculation now includes sex stratification. See “toxrefdb_2_0_recalc_pod.csv” for updated v2.0 POD values. The lowest loael/lel and highest noael/nel mg/kg/day value are selected for all corresponding studies, regardless of effect profile.

#-----------------------------------------------------------------------------------#
# v2.1 pods: lowest loael/lel, highest noael/nel
#-----------------------------------------------------------------------------------#

v2.1.pod <- dbGetQuery(con, "SELECT pod.*, study.study_type
                       FROM (prod_toxrefdb_2_1.pod
                       INNER JOIN prod_toxrefdb_2_1.study ON pod.study_id = study.study_id);") %>%
  data.table()

# extract study-level pod extremes
v2.1.pod.L <- v2.1.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(study_id, pod_type) %>%
   slice(which.min(mg_kg_day_value))
v2.1.pod.N <- v2.1.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(study_id, pod_type) %>%
   slice(which.max(mg_kg_day_value))
v2.1.pod.study <-rbind(v2.1.pod.N, v2.1.pod.L) %>%
  arrange(study_id, pod_type)

# extract chemical-level pod extremes
v2.1.pod.L <- v2.1.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(chemical_id, pod_type) %>%
   slice(which.min(mg_kg_day_value))
v2.1.pod.N <- v2.1.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(chemical_id, pod_type) %>%
   slice(which.max(mg_kg_day_value))
v2.1.pod.chem <-rbind(v2.1.pod.N, v2.1.pod.L) %>%
  arrange(chemical_id, pod_type)

#-----------------------------------------------------------------------------------#
# v2.0 pods: lowest loael/lel, highest noael/nel
#-----------------------------------------------------------------------------------#
v2.0.pod <- merge(read.csv("L:/Lab/NCCT_Toxrefdb/toxrefdb_2_1_public_release/toxrefdb_2_0_recalc_pod.csv"),dbGetQuery(con, "SELECT study.study_id, study.study_type FROM prod_toxrefdb_2_0.study;") , by="study_id") %>%
  data.table()

# extract study-level pod extremes
v2.0.pod.L <- v2.0.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(study_id, pod_type) %>%
   slice(which.min(mg_kg_day_value))
v2.0.pod.N <- v2.0.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(study_id, pod_type) %>%
   slice(which.max(mg_kg_day_value))
v2.0.pod.study <-rbind(v2.0.pod.N, v2.0.pod.L) %>%
  arrange(study_id, pod_type) %>%
  filter(study_id %in% v2.1.pod.study$study_id)

# extract chemical-level pod extremes
v2.0.pod.L <- v2.0.pod %>%
    filter(pod_type %in% c("loael", "lel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(qualifier, mg_kg_day_value)%>%
    group_by(chemical_id, pod_type) %>%
   slice(which.min(mg_kg_day_value))
v2.0.pod.N <- v2.0.pod %>%
    filter(pod_type %in% c("noael", "nel")) %>%
  filter(!is.na(mg_kg_day_value)) %>%
    arrange(desc(qualifier), mg_kg_day_value)%>%
    group_by(chemical_id, pod_type) %>%
   slice(which.max(mg_kg_day_value))
v2.0.pod.chem <-rbind(v2.0.pod.N, v2.0.pod.L) %>%
  arrange(chemical_id, pod_type) %>%
  filter(chemical_id %in% v2.1.pod.chem$chemical_id)  

#-----------------------------------------------------------------------------------#
# POD extremes data manipulation for version comparison
#-----------------------------------------------------------------------------------#

change.pod.study <- inner_join(v2.0.pod.study, v2.1.pod.study, suffix = c("_v2.0", "_v2.1"), by = c("study_id", "study_type", "pod_type")) %>%
  mutate(change = ifelse(mg_kg_day_value_v2.0 != mg_kg_day_value_v2.1, 1, 0)) %>%
  select(study_id, study_type, pod_type, mg_kg_day_value_v2.0, mg_kg_day_value_v2.1, change)

change.pod.chem <- inner_join(v2.0.pod.chem, v2.1.pod.chem, suffix = c("_v2.0", "_v2.1"), by = c("chemical_id", "pod_type")) %>%
  mutate(change = ifelse(mg_kg_day_value_v2.0 != mg_kg_day_value_v2.1, 1, 0)) %>%
  select(chemical_id, pod_type, mg_kg_day_value_v2.0, mg_kg_day_value_v2.1, change)

counts.pod.study <- change.pod.study %>%
  group_by(study_type, pod_type, change) %>%
  summarize(n_change = n())

counts.pod.chem <- change.pod.chem %>%
  group_by(pod_type, change) %>% 
  summarize(n_change = n())

podMetric.study.all <- change.pod.study %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change)) %>%
  select(study_id, pod_change)

podMetric.study.CHR <- change.pod.study %>%
  filter(study_type == "CHR") %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change))

podMetric.study.DEV <- change.pod.study %>%
  filter(study_type == "DEV") %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change))

podMetric.study.MGR <- change.pod.study %>%
  filter(study_type == "MGR") %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change))

podMetric.study.SAC <- change.pod.study %>%
  filter(study_type == "SAC") %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change))

podMetric.study.SUB <- change.pod.study %>%
  filter(study_type == "SUB") %>%
  group_by(study_id) %>%
  summarize(pod_change = sum(change))

podMetric.chem <- change.pod.chem %>%
  group_by(chemical_id) %>%
  summarize(pod_change = sum(change))

#-----------------------------------------------------------------------------------#
# examine magnitude of change
#-----------------------------------------------------------------------------------#
setDT(change.pod.study)
change.pod.study2 <- change.pod.study[change==1,]
change.pod.study2 <- change.pod.study2 %>%
  mutate(mg_kg_day_value_v2.1 = replace(mg_kg_day_value_v2.1, mg_kg_day_value_v2.1 == 0, 0.0000001)) %>%
  mutate(mg_kg_day_value_v2.0 = replace(mg_kg_day_value_v2.0, mg_kg_day_value_v2.0 == 0, 0.0000001))
change.pod.study2 <- change.pod.study2[, `:=` (magn.change = abs(mg_kg_day_value_v2.1 - mg_kg_day_value_v2.0), magn.change.log10 = abs(log10(mg_kg_day_value_v2.1) - log10(mg_kg_day_value_v2.0)))]

magn.PODchange.study.all <- change.pod.study2 %>%
  group_by(study_type, pod_type) %>%
  summarize(mean = mean(magn.change.log10), sd = sd(magn.change.log10), median = median(magn.change.log10), iqr = IQR(magn.change.log10), mad = mad(magn.change.log10))


setDT(change.pod.chem)
change.pod.chem2 <- change.pod.chem[change==1,]
change.pod.chem2 <- change.pod.chem2 %>%
  mutate(mg_kg_day_value_v2.1 = replace(mg_kg_day_value_v2.1, mg_kg_day_value_v2.1 == 0, 0.0000001)) %>%
  mutate(mg_kg_day_value_v2.0 = replace(mg_kg_day_value_v2.0, mg_kg_day_value_v2.0 == 0, 0.0000001))
change.pod.chem2 <- change.pod.chem2[, `:=` (magn.change = abs(mg_kg_day_value_v2.1 - mg_kg_day_value_v2.0), magn.change.log10 = abs(log10(mg_kg_day_value_v2.1) - log10(mg_kg_day_value_v2.0)))]

magn.PODchange.chem.all <- change.pod.chem2 %>%
  group_by(pod_type) %>%
  summarize(mean = mean(magn.change.log10), sd = sd(magn.change.log10), median = median(magn.change.log10), iqr = IQR(magn.change.log10), mad = mad(magn.change.log10))
#-----------------------------------------------------------------------------------#
# generate plots comparing study counts of pod changes between v2.0 and v2.1
#-----------------------------------------------------------------------------------#

# stacked style plot showing counts of studies with POD change and with no
# POD change, within study type/pod_type groups
change.study.stack <- ggplot(counts.pod.study, 
                             aes(x = pod_type, 
                                 y = n_change,
                                 label = n_change,
                                 fill = as.factor(change))) +
  geom_col(position = "stack") +
  geom_text(position = "stack", 
            vjust = -0.5, 
            size = 2.5) +
    scale_fill_manual(labels = c("POD unchanged", "POD changed"), 
                     values = c("aquamarine3", "darksalmon")) +
  labs(x = "POD type by study type", 
       y = "number of studies") +
  facet_grid(~study_type) +
  theme(axis.text.x = element_text(angle = 65, 
                                   vjust = 1, 
                                   hjust = 1), 
        aspect.ratio = 4/1 ) + 
  guides(fill = guide_legend(title = NULL))


change.study.stack

For all studies and study types, the number of changed values across POD types can be assessed. Options include: “No change”, or change in “One or more”, “Two or more”, “Three or more”, or “All four” POD types.

type <- c("all", "CHR", "DEV", "MGR", "SAC", "SUB")

No_change <- NULL
for(i in 1:length(type)){
  dft <- eval(parse(text = paste("podMetric.study", type[i], sep = ".")))
  No_change <- c(No_change, paste(round((nrow(dft[which(dft$pod_change == 0),])/n_distinct(dft$study_id))*100), "%", sep = ""))
}

One_or_more <- NULL
for(i in 1:length(type)){
  dft <- eval(parse(text = paste("podMetric.study", type[i], sep = ".")))
  One_or_more <- c(One_or_more, paste(round((nrow(dft[which(dft$pod_change >= 1),])/n_distinct(dft$study_id))*100), "%", sep = ""))
}

Two_or_more <- NULL
for(i in 1:length(type)){
  dft <- eval(parse(text = paste("podMetric.study", type[i], sep = ".")))
  Two_or_more <- c(Two_or_more, paste(round((nrow(dft[which(dft$pod_change >= 2),])/n_distinct(dft$study_id))*100), "%", sep = ""))
}

Three_or_more <- NULL
for(i in 1:length(type)){
  dft <- eval(parse(text = paste("podMetric.study", type[i], sep = ".")))
  Three_or_more <- c(Three_or_more, paste(round((nrow(dft[which(dft$pod_change >= 3),])/n_distinct(dft$study_id))*100), "%", sep = ""))
}

All_Four <- NULL
for(i in 1:length(type)){
  dft <- eval(parse(text = paste("podMetric.study", type[i], sep = ".")))
  All_Four <- c(All_Four, paste(round((nrow(dft[which(dft$pod_change == 4),])/n_distinct(dft$study_id))*100), "%", sep = ""))
}

Number_of_POD_types_changed <- c("All studies", "CHR", "DEV", "MGR", "SAC", "SUB")

Table2 <- data.frame(Number_of_POD_types_changed, No_change, One_or_more, Two_or_more, Three_or_more, All_Four)

colnames(Table2) <- c("Number of POD types changed", "No change", "One or more", "Two or more", "Three or more", "All four")

datatable(Table2,
          filter='top', 
          options=list(pageLength = 15, searching=FALSE, autoWidth=FALSE,  scrollX=TRUE, initComplete = JS(
            "function(settings, json) {", 
            "$('body').css({'font-family': 'Calibri'});", 
            "}")))

The distribution statistics for the magnitude of change within changed log10 POD values between v2.0 and v2.1, across all studies, are summarized in the below table.

grp.study <- rbind("CHR LEL", "CHR LOAEL", "CHR NEL", "CHR NOAEL", "DEV LEL", "DEV LOAEL", "DEV NEL", "DEV NOAEL", "MGR LEL", "MGR LOAEL", "MGR NEL", "MGR NOAEL", "SAC NOAEL*", "SUB LEL", "SUB LOAEL", "SUB NEL", "SUB NOAEL")
df.mag.study <- cbind(grp.study, round(magn.PODchange.study.all$median,2), round(magn.PODchange.study.all$iqr,2), round(magn.PODchange.study.all$mad,2), round(magn.PODchange.study.all$mean,2), round(magn.PODchange.study.all$sd,2))
df.mag.study <- rbind(df.mag.study, c("Overall", round(median(change.pod.study2$magn.change.log10), 2), round(IQR(change.pod.study2$magn.change.log10), 2), round(mad(change.pod.study2$magn.change.log10), 2), round(mean(change.pod.study2$magn.change.log10), 2), round(sd(change.pod.study2$magn.change.log10), 2)))

df.mag.study %>%
kable(col.names = c("group", "median", "IQR", "MAD", "mean", "sd"), caption = "Distribution statistics for magnitudes of change in log10 study-level PODs") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive")) %>%
  footnote(general = "SAC lel, SAC loael, and SAC nel are omitted since there are no changed PODs in these categories.",
           symbol = "No distribution spread since only 1 POD changed in this category.")
Distribution statistics for magnitudes of change in log10 study-level PODs
group median IQR MAD mean sd
CHR LEL 0.47 0.66 0.4 2.06 3.3
CHR LOAEL 0.3 0.43 0.26 1.49 3.02
CHR NEL 0.3 0.42 0.21 0.48 0.37
CHR NOAEL 0.16 0.26 0.19 0.27 0.3
DEV LEL 0.6 4.22 0.64 3.12 4.75
DEV LOAEL 0.56 4.58 0.76 3.27 5.15
DEV NEL 0.3 0.44 0.18 0.51 0.35
DEV NOAEL 0.3 0.1 0.06 0.36 0.18
MGR LEL 0.79 0.15 0.15 0.82 0.15
MGR LOAEL 4.74 8.59 6.45 4.93 5.28
MGR NEL 0.33 0.24 0.36 0.33 0.34
MGR NOAEL 0.35 0.63 0.41 0.44 0.46
SAC NOAEL* 0.3 0 0 0.3 NA
SUB LEL 3.81 7.51 5.59 3.82 4.36
SUB LOAEL 0.3 0.69 0.43 2.03 3.97
SUB NEL 0.54 0.45 0.35 0.52 0.33
SUB NOAEL 0.22 0.27 0.23 0.35 0.32
Overall 0.34 0.5 0.36 0.95 2.03
Note:
SAC lel, SAC loael, and SAC nel are omitted since there are no changed PODs in these categories.
* No distribution spread since only 1 POD changed in this category.
#-----------------------------------------------------------------------------------#
# plots showing distribution of study-level POD change magnitude by pod_type and study_type
#-----------------------------------------------------------------------------------#

#### box plots to show median and IQR
change.study.mag2 <- ggplot(data = change.pod.study2, 
                            aes(x = pod_type, y = magn.change.log10, fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_boxplot(outlier.shape = NA) +
  coord_cartesian(ylim = c(0, 10)) +
  scale_fill_manual("POD type", 
                    labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse", "darksalmon", "aquamarine3", "darkviolet")) +
  labs(x = element_blank(),
       y = "change in POD value (log10 mg/kg/day)") +
  facet_grid(~study_type) +
  theme_gray() + 
  theme(axis.text.x = element_blank(), 
        axis.text.y = element_text(size = 12),
        axis.ticks.x = element_blank()) + 
  guides(fill = guide_legend(title = NULL)) +
  ggtitle("Distributions of magnitudes of change in study-level PODs") +
  geom_text(
    data = change.pod.study2 %>% group_by(study_type, pod_type) %>% summarize(n = n()), 
    aes(y = -0.1, label = paste("n=", n, sep="")), size = 2 
  )

change.study.mag2

In the following visual, no density curves are produced for the SAC study type since there are no changed PODs for this study type except for a single NOAEL POD.

#### density plots to show distribution curve 
change.study.mag3 <- ggplot(data = change.pod.study2, 
                            aes(x = magn.change.log10, fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_density(alpha = 0.3) +
  scale_fill_manual("POD type", 
                    labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse", "darksalmon", "aquamarine3", "darkviolet")) +
  labs(y = "Density",
       x = "Change in POD value (log10 mg/kg/day)") +
  facet_rep_grid(rows = vars(study_type),
                  repeat.tick.labels = "bottom") +
  theme_gray() + 
  guides(fill = guide_legend(title = NULL)) +
  ggtitle("Distributions of magnitudes of change in study-level PODs")

change.study.mag3

#-----------------------------------------------------------------------------------#
# plots showing distribution of study-level POD change magnitude by study type and pod_type 
#-----------------------------------------------------------------------------------#

#### ridgeline density plots to show distribution curve 
change.study.mag4 <- ggplot(data = change.pod.study2, 
                            aes(x = magn.change.log10, y = pod_type)) +
  geom_density_ridges(aes(fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael"))), alpha = 0.3, show.legend = FALSE) +
  scale_fill_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse", "darksalmon", "aquamarine3", "darkviolet")) +
  labs(y = "Density",
       x = "Change in POD value (log10 mg/kg/day)") +
  theme_gray() + 
  ggtitle("Distributions of magnitudes of change in study level PODs overall")


change.study.mag4 + 
  geom_text(
    data = change.pod.study2 %>% group_by(pod_type) %>% summarize(n = n()), 
    aes(x = 4, label = paste("n=", n, sep="")), size = 3, vjust = -1.5 
  )

#-----------------------------------------------------------------------------------#
# plots showing distribution of study-level POD change magnitude by pod_type and study_type
#-----------------------------------------------------------------------------------#

#### simple plot of mean +/- sd 
change.study.mag1 <- ggplot(magn.PODchange.study.all, 
                             aes(x = pod_type, 
                                 y = mean,
                                 color = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_point() + 
  geom_errorbar(aes(ymin = mean - sd,
                    ymax = mean + sd)) +
  scale_color_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse", "darksalmon", "aquamarine3", "darkviolet")) +
  labs(x = element_blank(), 
       y = "mean change in POD value +- sd (log10 mg/kg/day)") +
  facet_grid(~study_type) +
  theme_gray() + 
  theme(axis.text.x = element_blank(), 
        axis.text.y = element_text(size = 12),
        axis.ticks.x = element_blank()) + 
  guides(fill = guide_legend(title = NULL))


change.study.mag1

3.4 Chemical-Level Changes in v2.0 to v2.1 PODs

There were 735 chemicals with calculated POD values in v2.1 to compare back to v2.0. The v2.0 pod table included calculated POD values for 1137 chemicals, each with 4 POD types. To do the comparison, the v2.1 POD calculation was rerun against v2.0 schema since calculation now includes sex stratification. See “toxrefdb_2_0_recalc_pod.csv” for updated v2.0 POD values. The lowest loael/lel and highest noael/nel mg/kg/day value was selected for all corresponding chemicals, regardless of effect profile.

#-----------------------------------------------------------------------------------#
# generate plots comparing chemical counts of pod changes between v2.0 and v2.1
#-----------------------------------------------------------------------------------#

# stacked style plot showing counts of chemicals with POD change and with no
# POD change, within pod_type groups
change.chem.stack <- ggplot(counts.pod.chem, 
                             aes(x = pod_type, 
                                 y = n_change,
                                 label = n_change,
                                 fill = as.factor(change))) +
  geom_col(position = "stack") +
  geom_text(position = "stack", 
            vjust = -0.5, 
            size = 2.5) + 
  scale_fill_manual(labels = c("POD unchanged", "POD changed"), 
                     values = c("aquamarine3", "darksalmon")) +
  labs(x = "POD type", 
       y = "number of chemicals") +
  theme(axis.text.x = element_text(angle = 65,
                                   vjust = 1, 
                                   hjust = 1),
        aspect.ratio = 3/4) + 
  guides(fill = guide_legend(title = NULL))


change.chem.stack

For all chemicals, the number of changed values across POD types can be assessed. Options include: “No change”, or change in “One or more”, “Two or more”, “Three or more”, or “All four” POD types.

No_change <- paste(round((nrow(podMetric.chem[which(podMetric.chem$pod_change == 0),])/n_distinct(podMetric.chem$chemical_id))*100), "%", sep = "")

One_or_more <- paste(round((nrow(podMetric.chem[which(podMetric.chem$pod_change >= 1),])/n_distinct(podMetric.chem$chemical_id))*100), "%", sep = "")

Two_or_more <- paste(round((nrow(podMetric.chem[which(podMetric.chem$pod_change >= 2),])/n_distinct(podMetric.chem$chemical_id))*100), "%", sep = "")

Three_or_more <- paste(round((nrow(podMetric.chem[which(podMetric.chem$pod_change >= 3),])/n_distinct(podMetric.chem$chemical_id))*100), "%", sep = "")

All_Four <- paste(round((nrow(podMetric.chem[which(podMetric.chem$pod_change == 4),])/n_distinct(podMetric.chem$chemical_id))*100), "%", sep = "")

Number_of_POD_types_changed <- "All chemicals"

Table2 <- data.frame(Number_of_POD_types_changed, No_change, One_or_more, Two_or_more, Three_or_more, All_Four)

colnames(Table2) <- c("Number of POD types changed", "No change", "One or more", "Two or more", "Three or more", "All four")

datatable(Table2,
          filter = 'top', 
          options = list(pageLength = 15, 
                         searching = FALSE, 
                         autoWidth = FALSE,  
                         scrollX = TRUE, 
                         initComplete = JS(
                           "function(settings, json) {", 
                          "$('body').css({'font-family': 'Calibri'});", 
                          "}")))

The distribution statistics for the magnitude of change within changed log10 POD values between v2.0 and v2.1, across all chemicals, are summarized in the below table.

grp.chem <- rbind("LEL", "LOAEL", "NEL", "NOAEL")
df.mag.chem <- cbind(grp.chem, round(magn.PODchange.chem.all$median,2), round(magn.PODchange.chem.all$iqr,2), round(magn.PODchange.chem.all$mad,2), round(magn.PODchange.chem.all$mean,2), round(magn.PODchange.chem.all$sd,2))
df.mag.chem <- rbind(df.mag.chem, c("Overall", round(median(change.pod.chem2$magn.change.log10), 2), round(IQR(change.pod.chem2$magn.change.log10), 2), round(mad(change.pod.chem2$magn.change.log10), 2), round(mean(change.pod.chem2$magn.change.log10), 2), round(sd(change.pod.chem2$magn.change.log10), 2)))

df.mag.chem %>%
kable(col.names = c("group", "median", "IQR", "MAD", "mean", "sd"), caption = "Distribution statistics for magnitudes of change in chemical-level PODs") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
Distribution statistics for magnitudes of change in chemical-level PODs
group median IQR MAD mean sd
LEL 0.45 0.73 0.4 1.11 1.85
LOAEL 0.42 0.68 0.37 1.28 2.35
NEL 0.28 0.35 0.27 0.37 0.34
NOAEL 0.24 0.32 0.21 0.36 0.41
Overall 0.3 0.5 0.26 0.77 1.55
#-----------------------------------------------------------------------------------#
# plots showing distribution of chemical-level POD change magnitude by pod_type 
#-----------------------------------------------------------------------------------#

#### box plots to show median and IQR since data is heavily skewed
change.chem.mag2 <- ggplot(data = change.pod.chem2, 
                            aes(x = pod_type, y = magn.change.log10,
                                fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_boxplot(outlier.shape = NA) +
  coord_cartesian(ylim = c(0, 2.5)) +
  scale_fill_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse",  "darksalmon", "aquamarine3","darkviolet")) +
  labs(x = element_blank(),
       y = "change in POD value (log10 mg/kg/day)") +
  theme_gray() + 
  theme(axis.text.x = element_blank(), 
        axis.text.y = element_text(size = 12),
        axis.ticks.x = element_blank()) + 
  guides(fill = guide_legend(title = NULL)) +
  ggtitle("Distributions of magnitudes of change in chemical-level PODs") +
  geom_text(
    data = change.pod.chem2 %>% group_by(pod_type) %>% summarize(n = n()), 
    aes(y = -0.05, label = paste("n=", n, sep="")) 
  )

change.chem.mag2

#-----------------------------------------------------------------------------------#
# plots showing distribution of chemical-level POD change magnitude by pod_type 
#-----------------------------------------------------------------------------------#

#### density plots to show distribution curve 
change.chem.mag3 <- ggplot(data = change.pod.chem2, 
                            aes(x = magn.change.log10, fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_density(alpha = 0.3) +
  scale_fill_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse",  "darksalmon", "aquamarine3","darkviolet")) +
  labs(y = "Density",
       x = "Change in POD value (log10 mg/kg/day)") +
  theme_gray() + 
  guides(fill = guide_legend(title = NULL)) +
  ggtitle("Distributions of magnitudes of change in chemical-level PODs")

change.chem.mag3

#-----------------------------------------------------------------------------------#
# plots showing distribution of chemical-level POD change magnitude by pod_type 
#-----------------------------------------------------------------------------------#

#### ridgeline density plots to show distribution curve 
change.chem.mag4 <- ggplot(data = change.pod.chem2, 
                            aes(x = magn.change.log10, y = pod_type)) +
  geom_density_ridges(aes(fill = factor(pod_type, levels = c("nel", "noael", "lel", "loael"))), alpha = 0.3, show.legend = FALSE) +
  scale_fill_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse", "darksalmon", "aquamarine3", "darkviolet")) +
  labs(y = "Density",
       x = "Change in POD value (log10 mg/kg/day)") +
  theme_gray() + 
  ggtitle("Distributions of magnitudes of change in chemical-level PODs")

change.chem.mag4 + 
  geom_text(
    data = change.pod.chem2 %>% group_by(pod_type) %>% summarize(n = n()), 
    aes(x = 4, label = paste("n=", n, sep="")), size = 3, vjust = -1.5 
  )

#-----------------------------------------------------------------------------------#
# plots showing distribution of chemical-level POD change magnitude by pod_type 
#-----------------------------------------------------------------------------------#

#### simple plot of mean +/- sd 
change.chem.mag1 <- ggplot(magn.PODchange.chem.all, 
                             aes(x = pod_type, 
                                 y = mean,
                                 color = factor(pod_type, levels = c("nel", "noael", "lel", "loael")))) +
  geom_point() + 
  geom_errorbar(aes(ymin = mean - sd,
                    ymax = mean + sd)) +
  scale_color_manual("POD type", 
                     labels = c("NEL", "NOAEL", "LEL", "LOAEL"), 
                     values = c("chartreuse",  "darksalmon", "aquamarine3","darkviolet"))+
  labs(x = element_blank(), 
       y = "mean change in POD value +- sd (log10 mg/kg/day)") +
  theme_gray() + 
  theme(axis.text.x = element_blank(), 
        axis.text.y = element_text(size = 12),
        axis.ticks.x = element_blank()) + 
  guides(fill = guide_legend(title = NULL))


change.chem.mag1

4 Reproducibility

4.1 File Export

This is code used to export the file with study-chemical level information accompanying the ToxRefDB v2.1 database release.

#write study-chem summary xlsx 
study_chem <- dbGetQuery(con, "SELECT *, study.* FROM prod_toxrefdb_2_1.chemical
                      INNER JOIN prod_toxrefdb_2_1.study on study.chemical_id=chemical.chemical_id;")

write.xlsx(study_chem, keepNA=TRUE, file='toxrefdb_2_1_study_chem_summary.xlsx')

4.2 Session

This is the R session environment and package information used to create this Rmarkdown release note.

print(sessionInfo())
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 19042)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=English_United States.1252 
## [2] LC_CTYPE=English_United States.1252   
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.1252    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] readr_2.1.2           patchwork_1.1.1       lemon_0.4.5          
##  [4] tidyr_1.2.1           kableExtra_1.3.4.9000 htmlTable_2.4.1      
##  [7] DT_0.23               plotly_4.10.0         ggridges_0.5.3       
## [10] ggplot2_3.3.6         dplyr_1.0.10          RMySQL_0.10.23       
## [13] DBI_1.1.3             data.table_1.14.4     openxlsx_4.2.5       
## [16] magrittr_2.0.3       
## 
## loaded via a namespace (and not attached):
##  [1] Rcpp_1.0.9        svglite_2.1.0     lattice_0.20-45   assertthat_0.2.1 
##  [5] digest_0.6.30     utf8_1.2.2        R6_2.5.1          plyr_1.8.7       
##  [9] backports_1.4.1   evaluate_0.18     highr_0.9         httr_1.4.4       
## [13] pillar_1.8.1      rlang_1.0.6       lazyeval_0.2.2    rstudioapi_0.14  
## [17] jquerylib_0.1.4   checkmate_2.1.0   rmarkdown_2.18    labeling_0.4.2   
## [21] webshot_0.5.4     stringr_1.4.1     htmlwidgets_1.5.4 munsell_0.5.0    
## [25] compiler_4.1.3    xfun_0.31         pkgconfig_2.0.3   systemfonts_1.0.4
## [29] htmltools_0.5.3   tidyselect_1.2.0  tibble_3.1.8      gridExtra_2.3    
## [33] fansi_1.0.3       viridisLite_0.4.1 tzdb_0.3.0        withr_2.5.0      
## [37] grid_4.1.3        jsonlite_1.8.3    gtable_0.3.1      lifecycle_1.0.3  
## [41] scales_1.2.1      zip_2.2.0         cli_3.3.0         stringi_1.7.8    
## [45] cachem_1.0.6      farver_2.1.1      xml2_1.3.3        bslib_0.4.1      
## [49] ellipsis_0.3.2    generics_0.1.3    vctrs_0.5.0       tools_4.1.3      
## [53] glue_1.6.2        purrr_0.3.4       crosstalk_1.2.0   hms_1.1.2        
## [57] fastmap_1.1.0     yaml_2.3.6        colorspace_2.0-3  rvest_1.0.3      
## [61] knitr_1.40        sass_0.4.2